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We study spatial pattern formation and energy localization in the dynamics 
of an anharmonic chain with quadratic and quartic intersite potential subject 
to an optical, sinusoidally oscillating field and a weak damping. The zone- 
boundary mode is stable and locked to the driving field below a critical forcing 
that we determine analytically using an approximate model which describes 
mode interactions. Above such a forcing, a standing modulated wave forms 
for driving frequencies below the band-edge, while a "multibreather" state 
develops at higher frequencies. Of the former, we give an explicit approximate 
analytical expression which compares well with numerical data. At higher 



forcing space-time chaotic patterns are observed. 
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I. INTRODUCTION 

The dynamics of classical anharmonic lattices displays a rich variety of features. Already 
the simplest models, as the one-dimensional chain of equal-mass oscillators with nearest- 
neighbors nonlinear forces, exhibit nontrivial solutions like anharmonic waves [TLB, discrete 
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solitons || and breathers |4j]. Due to its simplicity, one of the most widely studied examples 
in this class is the celebrated Fermi-Pasta-Ulam (FPU) model |J, where the interparticle 
potential is a simple fourth-order polynomial in the relative displacements. Many investi- 
gations have focused on the process of energy equipartition among phonons, after having 
fed the energy into long wavelength modes, whose instability leads to the generation of 
solitons ||. More recently the complementary case where the energy is placed into the 
highest frequency mode has been considered [|7|-f| . Above a certain energy threshold, which 
vanishes as 1/N, N being the number of oscillators, this mode becomes modulationally un- 
stable [11,12], leading to the growth of spatial modulations of the displacement field with 



a given finite wavelength. The subsequent evolution consists into the creation of local- 
ized structures (envelope solitons) which interact inelastically and coalesce in a few, large 
amplitude breathers [TIJ. These have however a finite lifetime and decay slowly until the 
asymptotic state of energy equipartition is attained. 

Up to now, processes of this kind have been studied mainly for Hamiltonian lattices, 
and a quite natural question is to ask how such phenomena are affected by both forcing 
and/or dissipative mechanisms that may arise from the presence of an external field and of 
the interaction with other degrees of freedom, respectively. These issues were addressed by 
Rossler and Page []13 |, who found that localized modes can indeed exist in a sinusoidally 



driven, undamped FPU chain. More recently, the same authors performed detailed studies 



of the optical generation of such excitations under the action of suitable impulsive fields [i ! 



Later on, investigations of models of driven-damped antiferromagnetic chains JT3] opened 
also the way to some experiments, whose outcomes are interpreted as manifestations of 
intrinsic localization [IB]. Theoretical studies of the parametrically driven, discrete, nonlin- 



ear Schrodinger equation [pl| and of coupled oscillator systems [|17| have also been recently 



undertaken and led to the discovery of periodic, quasiperiodic and even chaotic breathers. 

However, the relation between the phenomena observed in the Hamiltonian case with 
those appearing for forced-damped lattices has not yet been studied in detail. Indeed, 
we expect significative differences due to the creation of stationary states with nontrivial 



spatial structures, i.e. pattern formation [20,21]. In our context, having already studied the 



process of formation of stable localized structures arising from modulational instability in 
the conservative case || , we are strongly motivated to see how the presence of forcing and 
damping affects this process. To remain close to the Hamiltonian case, we restrict ourselves 
to the case of small damping. 

Various types of forcing are in principle possible, depending on the physical situation 
under study. However, a general requirement for localization is to excite band-edge modes. 
For Klein-Gordon lattices this is naturally realized using a spatially uniform driving field, 



which has been shown to induce interesting pattern formation phenomena [HJ . On the other 
hand, this forcing would be uneffective for FPU lattices, because, due to the symmetry 
of the Hamiltonian, the zero mode is decoupled (see below). Alternatively, since spatial 
localization appears from the instability of band-edge modes, we choose in this paper to 
drive the system at the zone boundary wavelength. Moreover, we consider the simplest case 



in which the driving field oscillates sinusoidally in time. 

In Section II we introduce the model and further comment on some of its features. After 
introducing a specific mode expansion of the displacement field, we obtain approximate 
equations of motion for such modes in the weak damping limit. Some details of the numerical 
simulations are presented as well as the quantities used to detect localization and to study 
spectral properties. 

Section III deals with the full characterization of the weak forcing solution where the 
zone boundary mode is locked to the driving field. By increasing the forcing parameter such 
a solution becomes unstable. Approximate values of the critical forcing and frequency are 
analytically computed and successfully compared with numerical results. 

Just beyond this instability more complex spatial patterns are attained by the system. 
For small enough frequencies a stable standing wave arises, which we describe analytically 
by solving for the fixed point of a suitable truncated mode expansion (Section IV). When 
the driving frequency exceeds a given resonant value a multibreather, spatially aperiodic 
state is instead observed in numerical simulations (Section V). 

Section VI is devoted to some conclusions and to a brief discussion of the transition to 
chaos, which is observed by increasing the forcing. 

II. THE FORCED AND DAMPED FPU MODEL 

The equations of motion of forced-damped FPU oscillator chain are 

u n = u n+1 + u n _i - 1u n + (u n+1 - u n f + (tt n _i - u n f - ju n + f cos(ut - nn), (1) 

where u n denotes the displacement of n-th oscillator with respect to its equilibrium posi- 
tion. Periodic boundary conditions, u n+ N = u n , are assumed, with N being the number 
of oscillators. Dimensionless units are used such that the masses, the linear and nonlinear 
force constants and the lattice spacing are taken equal to unity. The forcing and damping 
strengths are gauged by the parameters / and 7, respectively, and to is the driving frequency. 
As already mentioned in the Introduction, the choice of the forcing with the shortest wave- 
length is expected to favour the growth of localized excitations. Physically, the last term in 
Eq. ([]]) models the interaction of a uniform electric field applied to a chain of alternating 
opposite charges; it can in fact be written as /(— l) n cosut) [II3J. Let us stress again that the 
more widely studied case of uniform forcing |Tj| is not viable for models like (Ij) , because 
the zero mode is completely decoupled from the others as a consequence of the invariance 
of ([!]) under the transformation u n — > u n + const. 

In order to describe the forced oscillations of the system it is convenient to represent the 
displacement field in the form 
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- E [ a k e l ^ t+kn ^ + at k e -^'- te )] , (2) 



where a k are complex mode amplitudes and — n < k < n the corresponding wavenumber. 
Throughout the paper we will mainly focus our attention on the dynamics of the zone- 
boundary mode, which we refer to as 7r-mode for the sake of brevity. 

The equations of motion for the amplitudes a k s are obtained by substituting Eq. (g) into 
Eq. (|I|). Similarly to what is done for the undamped case M, a considerable simplification 
is achieved by neglecting higher-order harmonics that are produced by the cubic terms (the 
so-called rotating wave approximation). Moreover, in the limit of weak damping 7<w we 
find the following set of approximate equations 

- 2iua k - iu^a k = (u 2 k - u 2 )a k + 5 k ^f + 6 J2 ^Ji^i^aJ+^-fe, ( 3 ) 

where io\ = 2(1 — cosh), 5 k ^ is equal to one for k = tt and zero otherwise and 

Gg uq2 = -[1 + cos(gi + q 2 ) + cos(k - q 2 ) + cos(k - q x ) - cos k - cosgi 

- cos q 2 - cos(k — gi — q 2 )\ . 

Eqs. (|3|) correspond to the positive frequency projection in the base given in Eq. (|), the 
negative frequency ones being obtained by replacing k —>■ — k and taking the complex con- 
jugate. 

In the following, we aim at comparing the analytical predictions that can be drawn from 
the set of approximate equations (0) with the direct numerical simulations of model (^). 
We have integrated the equations of motion (|]) by means of a fourth-order Runge-Kutta 
algorithm with a step ranging between 10~ 2 and 10 -3 . We have always chosen an initial 
condition with all oscillators in their equilibrium position u n (0) = and small (of the order 
of 10~ 5 ) random Gaussian distributed initial velocities tt n (0). Several series of simulations 
have been performed for different values of the parameters uj and / and for fixed 7 = 0.1. 
The latter choice guarantees that the condition 7 <^ u holds in the resonance regions close to 
the band-edge, i.e. \u\ ~ c^v = 2 which is of main interest here. Furthermore, the resulting 
time-scales are reasonably short to allow a real-time analysis for chains as long as N = 512. 
Nonetheless, some of the results reported below have been checked also for another series of 
simulations performed with 7 = 0.01. 

We have monitored the energy density along the chain 
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h n = ~U n + - [(Un+i ~ U n ) + (u n - U„_i) J + - [(u n+1 - U n ) + {u n - U n _i) J , (4) 

as well as the spectrum of mode energies 

e k = \U k \ 2 + u 2 \U k \ 2 (5) 

where 



1 N /N 

U * = -nr T E u n ^ = V ( afe6 ^ + a -^ tUJt ) ^ 

VjV n =l Z 



is the amplitude of the k-th Fourier mode, which can be efficiently computed using a standard 
FFT routine. In the limit in which the a^'s are slowly varying with respect to the forcing 
or reach a stationary value, one finds that £& = Nuj 2 \a,k\ 2 is constant in time. On the other 
hand, when monitoring quantities like h n and/or u n it is convenient to observe the system 
at time intervals that are integer multiples of the driving period. 

III. MODULATIONAL INSTABILITY OF THE vr-MODE 

Let us begin by considering solutions where only the n mode is excited. Such solutions 
are numerically observed to exist and be stable below a certain critical forcing f cr that 
depends on both uj and 7. Indeed, in the typical simulation described at the end of the 
previous Section, all modes with k ^ tt damp out on a time-scale set by the value of 7 
while \a n \ rapidly grows and finally appoaches a constant value. Such asymptotic amplitude 
a n is obtained by solving for the stationary solution of Eqs. (H), which corresponds to an 
oscillation locked to the driving field with a constant phase lag. We get 
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a; 2 — 4 — 12|a 7r | 2 — i^a; 
and writing a^ = \a n \ exp(i9 n ) one gets 



(7) 



9 « = &tan U 2 "4-12KpJ ' (8) 

where the squared modulus z = {a^ 2 is the solution the cubic equation 

144z 3 - 24(4 - uj 2 )z 2 + [(4 - uj 2 ) 2 + 7 V] z = f . (9) 

One can easily ascertain that the latter admits a single real root only for \u\ < uj* where 
uj* ~ 2 + y / 37/2, while three distinct roots may otherwise exist (see the discussion in the 
following) . 

Before treating in detail the differences between the two cases, let us address the question 
of stability. This is accomplished by solving the set of equations obtained linearizing Eqs. @ 
around the stationary solutions (i.e. neglecting all interaction terms which are nonlinear in 
the afcs) 

— 2iujdk — iu>jak = (u>l — to 2 )ak + S^la^ a~^ k , (10) 

where 

^(l + eMV 2 (11) 

is the frequency of the fc-th mode shifted by the interaction with the 7r-mode. As usual, Eq. 
(p!0|), together its complex conjugate, are solved looking for solutions of the form exp(z/fci). 



The relevant branch of the eigenvalue spectrum reads 



*k = ~\ + ^y/*SkK\*-(Z 2 k-u 2 ) 2 ■ (12) 

The growth rate Re{vk} is maximal when the square root in the above expression attains its 
maximum value, i.e. when the resonance condition uo = u)k t holds |22|]. The latter, together 
with definition (|TT|) , fixes the value of the wavenumber fc* of the most unstable mode as 

C ° SK = 1 -2(l + 6\a^) ■ (13) 

In the case in which a single real root of (^)) exists (i.e. \u>\ < a;*), the threshold for 
modulational instability can be computed explicitly. Indeed, by letting v k = in Eq. (|12"D 
and using formula (|i~3"D, one gets 

'^ = 3(uj - 2 7 ) ■ (M) 

Finally, from formula (0) one derives the value of the critical forcing 
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+ 7 2 w 2 } . (15) 



In the case where three solutions exist (\u>\ > a;*), the stability properties must be 
considered separately for each of them. Let us discuss this issue with reference to the case 
uj = 2.4, illustrated by the graph in Fig. Q. Here the three solutions, which we label A, B 
and C, coexist in the range / + < / < /_. The values of f± can be computed from Eq. @, 
since they correspond to the amplitudes 
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Hence, from formula (|^) we obtain 



(16) 



f ± = KW(u 2 - 4 - 12|a±| 2 ) 2 + 7 2 ^ 2 • (17) 

Looking at the spectrum (|T2"D, it turns out that the largest amplitude solution C is always 
modulationally unstable. The intermediate amplitude solution B is modulationally unstable 
above {a^l 2 = (u 2 — 4)/24 (corresponding to / = 0.235 in Fig. P while it has the maximal 
growth rate at ±7r below this point. The smaller amplitude solution A is instead found to be 
always stable. Since we choose to run the dynamics starting always with almost zero mode 
amplitudes, it would be natural to conclude that the resulting trajectory converges to the 
solution A up to its existence boundary and hence that f cr = /_. The numerical simulations 
show that this conclusion is not actually correct. Indeed, the system approaches the stable 
solution A only up to a critical value of the forcing /^* which is definitely smaller than /_ 
(see the diamonds in Fig. [1|). Beyond such a value, the 7r-mode undergoes modulational 
instability similar to the previous case. 



An interpretation of this phenomenon can be given in the following way. Let us consider 
the equations of motion for the "internal" dynamics of the 7r-mode that can be derived from 
Eq. (g), 

— 2icud n — iwya-n = (4 — u 2 )a n + f + 12a 7r |a 7r | 2 . (18) 

This nonlinear equation provides a good approximation of the dynamics as long as all the 
other modes are not significantly excited. Its numerical solution indicates that for / = 
/cr*(^) the initial condition a n (0) = 0,6^(0) = exits the basin of attraction of the fixed 
point corresponding to the solution A. Such a value corresponds pretty well to the actual 
critical forcing numerically determined for the full system (see the solid vertical line at /*"* in 
Fig. [I]). In other words, at /*"* the dynamics leaves the lowest amplitude solution A because 
of this internal instability and "jumps" into the modulationally unstable region. Since a 
fixed point will not be subsequently approached, we do not expect that the corresponding 
spectrum of growth rates will be described by formula QT2"!). Nonetheless, a reasonable 
qualitative agreement between the two is observed, being both characterized by a sharp 
maximum around the most unstable mode with some broadened band around it. 

Finally, the results described above are summarized in Fig. |2], where we show the control 
parameter space (f,u>). The dashed line for \u\ < u* is the theoretical expression ( fl5|) for 
the modulational instability valid when only one solution is present. At cj* two full lines 
start, given by formula (|i~7|), which bound the region where three solutions occur. Inside this 
region, the dotted line /*"* is the one obtained numerically by looking at the internal 7r-mode 
instability, as we have described above. The numerical data (full triangles) were obtained by 
looking at the incipient modulational instability of the full system. An excellent agreement 
(within some percent) with the theoretical results is observed in the considered frequency 
range. We have also checked that the critical wavenumber k* is accurately predicted by 
Eq. (0). We thus conclude that both the theory developed for the case of one solution and 
the approximate description of the internal instability when three solution are present are 
basically correct. 

IV. STANDING NONLINEAR WAVES 

Let us now describe the states forming just above the threshold f cr , after the development 
of the modulational instability. Two different behaviors appear depending whether the 
driving frequency lies below or above a resonance frequency, which for small 7 is very close 
to the upper band edge uj ~ 2. In this Section we discuss the first case. As expected from 
the spectrum fll2"D, the development of the instability leads to a fast growth of the modes 
belonging to the unstable band around k* (and also its harmonics, see below). Afterwards, 
the main band shrinks until the system saturates to the asymptotic state. The result is 
basically a modulated standing wave locked in time with the forcing field (see the example 
illustrated in Fig. |^) and the wavenumber of the modulation is indeed very close to the 
expected value ti — k* with k* given by Eq. dl3|) . Furthermore, this state appears to be 



stable, at least on the time scale considered in the simulations. For instance, the wave 
profile shown in Fig. [3] remains unaltered up to 2 10 5 time units, i.e. for more than 6 10 4 
driving periods. 

An approximate theoretical analysis can explain the formation of this pattern. Indeed, 
in view of the above results, it is reasonable to look for a simplified description that neglects 
all the modes but the 7r-mode and the most unstable one with wavenumber fc*. Under such 
assumptions and taking into account the resonance condition we obtain from Eqs. (|3|) the 
following coupled equations for a n and a kt : 

—2iud w — iwya^ = (4 — uo ja^ + f + 12a w \a w \ + 6co kt a^a k:IL a- kt + 12u k ^a ir \ak ti \ ■ 

9 
—2iuja kt — itujcbk, = 3w^oJoi fcj> + -j-ci^OkJajfeJ . (19) 



The first equation is nothing but the modified version of Eq. ( |T8"D with the interaction 
terms between it and k* modes included. The stationary solutions are determined letting 
a k, = \ a k t \ ex p(^fc») an d first solving the second of Eqs. (|19|) 

4 / la |2 

o - 2~V W 4 ~~ \ a -* \t > sin2(0 w - 9 ±k J = - — 
duoi v a 



a± k f - ^-Vlo^rl 4 - K\% , sin2%-9 ±kt ) = -^§ . (20) 



Substituting the latter in the first of Eqs. ( fL9|) we get the stationary value of |a,r| above 
threshold 
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f, (21) 



70; + 8y la^l 4 — \a 

Solving this equation, we can thus get \a± k *\ an d 9± ki from Eq. (|20|). A fourth equation, 
which we do not explicitly display here for the sake of brevity, allows to compute 6 n as well. 
The stationary solution obtained above corresponds to a displacement field which can 
be derived using expansion (g) 

u n = l a 7r|( — l) n cos(u;t + 6 n ) + 2\a kt \ cos(/c*n) cos(uot + 9 kt ) , (22) 

which is precisely a standing modulated wave. These analytical results are also in reasonable 
quantitative agreement with numerical data. For example, in the case of Fig. ^] (uo = 1.8, 
/ = 0.150 and 7 = 0.1) we get 

(23) 
(24) 

the relative deviations being as expected of the order j/00. 

The above expressions are simplified close to threshold and sufficiently far from resonance. 
Indeed in this case 9 n is vanishingly small (see Eq. (KJ)) and therefore one gets, from formulas 
(POD, 9 kt ~ 7r/4 and the approximate solution 
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a n \ 2 


= 0.02085 


l«fcj 2 


= 4.88 10 -4 (numerical) 


a^\ 2 


= 0.02085 


K,| 2 


= 4.2710 -4 (theoretical) 



7T 

u n ~ |o^|(— l) n cos(ujt) + 2 1 a*. J cos(fc*n) cos(cjt + — ) . (25) 

This solution compares well with numerical data. Looking at the pattern at times which are 
integer multiples of 2ti/uj, again for the example of Fig. [^ with the values in formula (p4]). 

u n ~ 0.1437(-l) n + 0.0292 cos(2.05n) (theoretical) 

u n = 0.142(-l) n + 0.0359 cos(2.025n) (numerical) , (26) 

which shows a good agreement. 

Before concluding this Section, let us discuss the issue of higher-order corrections to the 
solution (^). Indeed, besides the main component at k*, the nonlinear terms induce the 
presence of several (exponentially small) harmonics whose wavenumber can be expressed as 
(recall that — ix < k < n) 

k n = k*n + (n — l)n, n = 2, 3.... (27) 

Their amplitudes can be computed perturbatively from the stationary solution of Eq. @. 
For instance, the first harmonic (n = 2) is evaluated as a function a^ and a kt . We give here 
for completeness the explicit expressions of the first and second harmonics. 

18 cos K (cos K. - 1) 2 

l a fc 2 | — 5 o l a 7r||a/cj , (2o) 

3 1 1 — 6 cos 2 fc* + 2 cos 3 fc* I o 

a fc3 = r^ ^ \a k ,\ ■ (29) 

Higher-order harmonics (n > 2) can be thus computed recursively in a similar way obtaining 
the general result 

a kn = g _^ 2 51 G ^,P„ a <?n a Pn a -r n) 9n +Pn + T n = k n . (30) 

kn q n ,Pn,r n <k n 

Notice that an infinite number of harmonics is expected here. This shows that the waves 
we are dealing with are more general than those previously reported in the literature for 
damped-driven Klein-Gordon lattices [BJ. In this latter case, the modulation has a finite 



number of harmonics due to the peculiar mutual relationship among modes 0, 7r/2 and tc 
that allows for solution where no other modes are excited. 

V. MULTIBREATHER STATES 

A different scenario is observed for driving frequency above resonance, which we briefly 
describe here with reference to the case uj = 2.4, illustrated in Fig. ||. The instability 
described at end of Section III produces, on a relatively short time scale (t ~ 10 2 ), a 
disordered assembly of sharply localized structures in a similar way to what is observed for 
the undriven case |||Kj. On longer time scales (t ~ 10 3 ), a further stage follows in which 

9 



the localized peaks arrange themselves until they eventually reach an asymptotic state (see 
the upper panel of Fig. |j). This is a sort of "multibreather" state, i.e. an array of unevenly 
spaced breathers. Remarkably, such a complex solution appears to be stable as the previously 
illustrated nonlinear wave. For instance, the pattern in Fig. |] has been observed to persist 
up to a time 2.0 10 5 , i.e. more than 7 10 4 driving periods. 

The corresponding displacement field in Fig. || reveals that each breather is pretty similar 



to the "even-parity" localized modes found in the Hamiltonian case ||12| , the main difference 
being the presence of a 7r-mode background induced by the field. No "odd-parity" mode is 
generated by the above mechanism in all the examined cases. Actually, a closer inspection 
reveals that two distinct solutions of slightly different amplitude are present. Moreover, 
apart from a small lag induced by the damping, the background oscillates always out of 
phase with respect to the field as it can be ascertained by comparing the dotted line in 
Fig. |5| with the corresponding displacement pattern. On the contrary, the breathers are 
always in-phase with the field. Accordingly, depending on the site on which they sit, they 
can have equal or opposite relative phases. Although localized solutions oscillating out of 
phase with the field are known to exist in similar models [|13| |, they are not detected in the 
present context. 

The inhomogeneous distribution of vibrational energy is distinctly reflected in the mode 
spectrum (see the lower panel of Fig. |j), which displays a localized structure and a band 
broadening as a consequence of the irregular spacing among the breathers. 

To give a more firm basis to our numerical observations, it is useful to briefly point out 
some conclusions that can be drawn by a suitable continuum approximation of the FPU 
model. In analogy with the approach followed for the hamiltonian case |TD| , one can in fact 
write the displacement field as u n = (—l) n Re{ip(x,t)exp(iut)} where ip is an envelope 
function which is assumed to be slowly varying in space and time on scales of the order of 
the interatomic spacing and the driving period, respectively. A standard calculation leads 
to the driven-damped nonlinear Schrodinger equation 

2iui> + (4 - u 2 + iuj)ijj + </>,, + 12?/#| 2 = / , (31) 



which in the spatially uniform case reduces to Eq. (|18| ) for the 7r-mode amplitude. Notice 
however that this type of description makes sense only when a relatively narrow packet 
of modes close to ±7r is excited and is thus less general than the one based on Eqs. (^]). 
For a suitable choice of parameters, Eq. (|3~TD admits two soliton solutions of different am- 
plitudes |23] as well as stable multisoliton complexes arising from bifurcation of one of 



them [p4| . Thus, the localized states of Fig. [^ could be related to such solutions, at least 



to the extent in which our lattice model can be approximated by a continuum equation like 



(|3T|). Although a more quantitative comparison would be desirable, this is a solid argument 



in support of the existence of stable multibreather complexes. 
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VI. CONCLUSIONS AND PERSPECTIVES 

We have confirmed that modulational instability of zone boundary modes is a relevant 
mechanism for the generation of nontrivial spatial structures in discrete anharmonic lattices. 
Our results for the externally driven case complement previous studies on Hamiltonian mod- 
els and further show that the interplay of (almost) resonant forcing and damping can stabilize 
such structures. Our approximate analytical framework has allowed also to derive the stabil- 
ity chart of the zone-boundary mode. For forcing frequencies below resonance a modulated 
wave is formed after instability, whose approximate analytical expression we have derived in 
terms of mode amplitudes. For frequencies above resonance a "multibreather" state arises 
of which we have given a phenomenological characterization. There are other problems that 
could be attacked within this approach, like the interesting issue of destabilization of the 
modulated wave. This could also contribute to a better understanding of localization mech- 
anisms. Furthermore, this points out the possibility of generating long-living and complex 
energy distributions in space for real ionic crystals under the action of optical fields. 

Let us also briefly comment on the fate of the two typical spatial patterns described 
above. As expected, increasing the forcing leads to their destabilization and eventually 
to the transition to chaos (see the stars in Fig. ||] as well as Fig. ^). This bears a strong 
analogy with the parametrically driven-damped nonlinear Schrodinger equation where stable 
patterns precede the onset of chaotic behavior p5 |. This is a further indication that a 



considerable insight on the dynamics of the discrete model may be achieved from Eq. (pip . 
Besides of this, it has been suggested on quite general grounds that systems obeying the 
symmetry «„—>«„ + const, should show an abrupt transition to a soft-turbulent state 
without displaying any intermediate spatial pattern |j26f . This is in apparent contradiction 
with the existence of stable patterns found in our numerics. On the other hand, performing 
simulations at 7 = the transition to chaos occurs as soon as / 7^ 0, since the modulational 
pattern does not attain any regular asymptotic state. We can therefore conjecture that the 
presence of friction alters the nature of the transition in a way that remains to be understood. 
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FIG. 1. Squared amplitude of the multiple fixed-point solutions vs. the forcing / for uj = 2.4 
and 7 = 0.1. The solid vertical line represents the critical value of the forcing f%r = 0.224. Dashed 
lines are the boundaries for the range /+</</_ where three solutions exist. The letters A,B,C 
denote the three solutions from the smallest to the highest amplitude one. 
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FIG. 2. Control parameter plane (u>, /) for 7 = 0.1. The dashed line for uj < lo* is the critical 
forcing f CT given by formula (|lq) at which modulational instability occurs. The region for uj > to* 
where three solutions exist is bounded by the solid lines f±. The dotted line within corresponds 
to the instability threshold /*"* evaluated numerically from Eq. (|18|). The full triangles are the 
numerical estimates of f cr . Open circles left (resp. right) of the to = u* line denote points where 
standing waves (resp. breathers) occur. The stars are some parameter values for which chaos is 
detected. 
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FIG. 3. Pattern generated after the modulational instability for u> = 1.8, / = 0.150 
(/cr = 0.148). This pattern stabilizes at t — 10 4 . In the upper panel we show the energy density in 
a part of a chain of 512 particles, in the lower the corresponding mode energy spectrum in lin-log 
scale. 
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energy density, in the lower the mode energy spectrum. 
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FIG. 6. Chaotic state obtained for uj = 1.8, / = 0.27. In the upper panel we show a snapshot 
of the energy density along the chain at t = 5200 and in the lower the corresponding mode energy 
spectrum in lin-log scale. 



